Using Rdata file as of 13-04-22
mood_data = as_tibble(get(load('my_data_Questionnaire.Rdata'))) %>%
filter(questionnaire == 'Checklist')
# # filter for participants with 4 sessions
# completed = c()
# for (i in 1:length(unique(mood_data$participant))) {
# ppt = unique(mood_data$participant)[i]
# by.participant = filter(mood_data, participant == ppt)
# N.sessions = length(unique(by.participant$session))
# if (N.sessions == 4) {
# completed = append(completed, ppt)
# }
# }
# setdiff(completed, mylist)
# this is the list of valid ppts to include in analysis
# from gender.xlsx in lagringshotell (-E019), n=61
completed = c(1,2,3,4,5,8,10,12,14,15,17,19,21,22,23,24,26,29,30,32,33,34,36,37,
38,39,43,44,45,46,47,51,53,55,58,59,60,61,62,63,64,66,71,73,75,77,
78,80,81,84,86,87,90,92,97,99,100,103,107,112,120,123)
mood_data_c = mood_data %>%
filter(participant %in% completed) %>%
mutate(task_order = ifelse(stress_sessions == '1_3', 1, 0)) %>% # 1 = stress first, 0 = control first
select(-c(status, pilot, comment, questionnaire, item_order, RT, date, gender,
drug_sessions, stress_stim, pleas_stim, stress_sessions)) %>%
mutate(item_name = recode(item_name, "conscious" = "self_conscious",
"desperate" = "distressed")) %>%
dplyr::rename(task_type = stress_type)
# # particpant 29 session 1 didn't get read into Guido's dataset, so I'm manually adding it here
# E029_1 = read.csv('E029_1_2022_Jan_10_0953_questionnaire_custom.csv') %>%
# as_tibble() %>%
# filter(questionnaire == 'Checklist') %>%
# mutate(Drug = 'placebo', Stress = 'stress', task_type = 'tsst', task_order = 1) %>%
# select(-c(questionnaire, item_order, item_text, RT)) %>%
# mutate(participant = as.numeric(recode(participant, 'E028' = '29'))) %>%
# mutate(stage = as.numeric(recode(stage, 'Q1'='1', 'Q2'='2', 'Q3'='3', 'Q4'='4','Q5'='5',
# 'Q6'='6', 'Q7'='7', 'Q8'='8', 'Q9'='9', 'Q10'='10', 'Q11'='11'))) %>%
# mutate(item_name = recode(item_name, "conscious" = "self_conscious",
# "desperate" = "distressed"))
#
# mood_data_c = bind_rows(mood_data_c, E029_1)
mood_data_c$participant = as.factor(mood_data_c$participant)
mood_data_c$session = as.factor(mood_data_c$session)
mood_data_c$stage = as.factor(mood_data_c$stage)
#mood_data_c$item_order = as.factor(mood_data_c$item_order)
mood_data_c$item_name = as.factor(mood_data_c$item_name)
#mood_data_c$gender = as.factor(mood_data_c$gender)
mood_data_c$Drug = as.factor(mood_data_c$Drug)
mood_data_c$Stress = as.factor(mood_data_c$Stress)
mood_data_c$task_type = as.factor(mood_data_c$task_type)
mood_data_c$task_order = as.factor(mood_data_c$task_order)
#mood_data_c$stress_sessions = as.factor(mood_data_c$stress_sessions)
#mood_data_c$drug_sessions = as.factor(mood_data_c$drug_sessions)
#mood_data_c$stress_stim = as.factor(mood_data_c$stress_stim)
#mood_data_c$pleas_stim = as.factor(mood_data_c$pleas_stim)
# mood_data_c = mood_data_c %>%
# arrange(participant, session, stage) %>%
# ddply(.(participant, session, item_name), transform, bsln_corr = response - response[3])
mood_data_c_placebo = mood_data_c %>%
filter(Drug == "placebo")
Make some plots
mood_data_c %>%
ggplot(aes(x = stage, y = response, color = Stress, group = participant:session)) +
geom_line() +
facet_wrap(~item_name) +
theme(legend.position = "top") +
geom_vline(xintercept = 3.5, col = "magenta", lty = 2) +
geom_vline(xintercept = 4.5, col = "blue", lty = 2) +
geom_vline(xintercept = c(5.5,7.5), col = "black", lty = 2)+
labs(title = "All sessions")
mood_data_c_placebo %>%
ggplot(aes(x = stage, y = response, color = Stress, group = participant:session)) +
geom_line() +
facet_wrap(~item_name) +
theme(legend.position = "top") +
geom_vline(xintercept = 3.5, col = "magenta", lty = 2) +
geom_vline(xintercept = 4.5, col = "blue", lty = 2) +
geom_vline(xintercept = c(5.5,7.5), col = "black", lty = 2)+
labs(title = "Placebo only")
# raw_means = mood_data_c %>%
# summarySEwithin(measurevar = 'response', withinvars = c('Stress', 'stage', 'item_name'), idvar = 'participant') %>%
# mutate(lower = response-2*se, upper = response+2*se)
raw_means = mood_data_c %>%
summarySEwithin2(measurevar = 'response', withinvars = c('Stress', 'stage', 'item_name'), idvar = 'participant') %>%
mutate(lower = response-ci, upper = response+ci)
## Warning in qt(conf.interval/2 + 0.5, unlist(datac$N) - 1): NaNs produced
## Warning in qt(conf.interval/2 + 0.5, unlist(datac$N) - 1): NaNs produced
raw_means %>%
ggplot(aes(x = stage, y = response, color = Stress, group = Stress)) +
geom_line() +
geom_ribbon(aes(ymin = lower, ymax = upper, color = Stress, fill = Stress), alpha = .2) +
facet_wrap(~item_name)+
labs(title = "All sessions")
raw_means_placebo = mood_data_c_placebo %>%
summarySEwithin2(measurevar = 'response', withinvars = c('Stress', 'stage', 'item_name'), idvar = 'participant') %>%
mutate(lower = response-ci, upper = response+ci)
## Warning in qt(conf.interval/2 + 0.5, unlist(datac$N) - 1): NaNs produced
## Warning in qt(conf.interval/2 + 0.5, unlist(datac$N) - 1): NaNs produced
raw_means_placebo %>%
ggplot(aes(x = stage, y = response, color = Stress, group = Stress)) +
geom_line() +
geom_ribbon(aes(ymin = lower, ymax = upper, color = Stress, fill = Stress), alpha = .2) +
facet_wrap(~item_name)+
labs(title = "Placebo only")
raw_means = mood_data_c %>%
summarySEwithin2(measurevar = 'response', withinvars = c('task_type', 'stage', 'item_name'), idvar = 'participant') %>%
mutate(lower = response-ci, upper = response+ci)
## Warning in qt(conf.interval/2 + 0.5, unlist(datac$N) - 1): NaNs produced
## Warning in qt(conf.interval/2 + 0.5, unlist(datac$N) - 1): NaNs produced
raw_means %>%
ggplot(aes(x = stage, y = response, color = task_type, group = task_type)) +
geom_line() +
geom_ribbon(aes(ymin = lower, ymax = upper, color = task_type, fill = task_type), alpha = .2) +
facet_wrap(~item_name)+
labs(title = "All sessions")
raw_means_placebo = mood_data_c_placebo %>%
summarySEwithin2(measurevar = 'response', withinvars = c('task_type', 'stage', 'item_name'), idvar = 'participant') %>%
mutate(lower = response-ci, upper = response+ci)
## Warning in qt(conf.interval/2 + 0.5, unlist(datac$N) - 1): NaNs produced
## Warning in qt(conf.interval/2 + 0.5, unlist(datac$N) - 1): NaNs produced
raw_means_placebo %>%
ggplot(aes(x = stage, y = response, color = task_type, group = task_type)) +
geom_line() +
geom_ribbon(aes(ymin = lower, ymax = upper, color = task_type, fill = task_type), alpha = .2) +
facet_wrap(~item_name)+
labs(title = "Placebo only")
These have a few missing values: participant 58, session 3, stage 9 participant 5, session 2, stage 11
mood_efa_withdrug = mood_data_c %>%
filter(item_name %!in% c("stomach_discomfort", "dry_mouth", "dizzy",
"nauseous", "high", "euphoric", "venflon_pain")) %>%
pivot_wider(names_from = item_name, values_from = response)
save(mood_efa_withdrug, file = "mood_efa_withdrug.Rdata")
mood_efa_nodrug = mood_data_c %>%
filter(item_name %!in% c("stomach_discomfort", "dry_mouth", "dizzy",
"nauseous", "high", "euphoric", "venflon_pain")) %>%
filter(Drug == 'placebo') %>%
pivot_wider(names_from = item_name, values_from = response)
save(mood_efa_nodrug, file = "mood_efa_nodrug.Rdata")
Extraction: Maximum likelihood Factor retain: parallel Rotation: Oblimin
jmv::efa(
data = mood_efa_nodrug,
vars = vars(bored, confident, distressed, dull, embarassed, good, happy,
heart_palpitations, indifferent, irritable, not_myself, relaxed,
safe, shaking, stressed, vulnerable, warm_face, anxious),
extraction = "ml",
hideLoadings = 0.0,
sortLoadings = TRUE,
screePlot = TRUE,
eigen = TRUE,
factorCor = TRUE,
factorSummary = TRUE,
modelFit = TRUE)
## Loading required namespace: GPArotation
##
## EXPLORATORY FACTOR ANALYSIS
##
## Factor Loadings
## -------------------------------------------------------------------------------------------------------------------
## 1 2 3 4 5 Uniqueness
## -------------------------------------------------------------------------------------------------------------------
## happy 0.947730229 -0.002839789 0.03186217 -0.010008481 -0.103344116 0.1441011
## good 0.880248990 -0.017255961 -0.01737104 -0.044592733 0.007230154 0.1964062
## confident 0.515606180 0.021851332 -0.03101663 0.116687837 0.367740950 0.4899410
## relaxed 0.493891642 -0.083842036 -0.10931276 0.140841342 0.198064631 0.5369914
## safe 0.489935964 -0.081674889 -0.03814392 -0.012103937 0.349386788 0.4895547
## distressed -0.037512388 0.837165324 -0.02355960 -0.012455925 0.117733891 0.3496694
## anxious -0.026918046 0.556211880 0.19190670 0.025816950 -0.089693050 0.4430262
## embarassed -0.037466082 0.541435265 0.13811042 -0.010166279 -0.214798747 0.4271340
## vulnerable -0.006539164 0.418706723 0.20111668 0.086167236 -0.330636258 0.4457756
## stressed -0.024289410 0.401826325 0.33850415 -0.092969900 -0.196639473 0.3778799
## irritable -0.091165560 0.382365381 0.08975423 0.169106225 0.115345828 0.7514648
## not_myself 0.017914057 0.275482822 0.22618986 0.053453244 0.141570170 0.8042805
## heart_palpitations -0.009598171 0.020001116 0.73755628 -0.086114450 -0.015097044 0.4116822
## warm_face 0.033170237 0.024888705 0.65860149 0.023250675 0.093470054 0.5732425
## shaking -0.036582665 0.167791649 0.64007186 -0.007495731 -0.057350108 0.3644579
## bored -0.058062399 0.080329488 -0.10516874 0.709371544 0.046611826 0.4676823
## indifferent 0.083839042 -0.030680056 -0.03318914 0.655874860 -0.191656258 0.5864595
## dull -0.117636934 -0.145578113 0.34318913 0.464949594 0.204042116 0.6383806
## -------------------------------------------------------------------------------------------------------------------
## Note. 'Maximum likelihood' extraction method was used in combination with a 'oblimin' rotation
##
##
## FACTOR STATISTICS
##
## Summary
## ----------------------------------------------------------
## Factor SS Loadings % of Variance Cumulative %
## ----------------------------------------------------------
## 1 2.6414143 14.674524 14.67452
## 2 2.4880064 13.822258 28.49678
## 3 2.2692575 12.606986 41.10377
## 4 1.2567800 6.982111 48.08588
## 5 0.8464121 4.702289 52.78817
## ----------------------------------------------------------
##
##
## Inter-Factor Correlations
## -------------------------------------------------------------------
## 1 2 3 4 5
## -------------------------------------------------------------------
## 1 — -0.3607216 -0.2653469 -0.05073137 0.2176539
## 2 — 0.7065712 -0.04512222 -0.3065451
## 3 — -0.04305492 -0.2301821
## 4 — 0.2060010
## 5 —
## -------------------------------------------------------------------
##
##
## MODEL FIT
##
## Model Fit Measures
## ----------------------------------------------------------------------------------------------------
## RMSEA Lower Upper TLI BIC <U+03C7>² df p
## ----------------------------------------------------------------------------------------------------
## 0.05831881 0.05288600 0.06393751 0.9307498 -118.7192 407.2379 73 < .0000001
## ----------------------------------------------------------------------------------------------------
##
##
## EIGENVALUES
##
## Initial Eigenvalues
## -------------------------
## Factor Eigenvalue
## -------------------------
## 1 5.42307764
## 2 1.43170852
## 3 0.93422306
## 4 0.27000197
## 5 0.14056087
## 6 0.09026108
## 7 0.04381429
## 8 -0.03579408
## 9 -0.09506274
## 10 -0.11636128
## 11 -0.18831231
## 12 -0.19091637
## 13 -0.21736653
## 14 -0.26051196
## 15 -0.34816400
## 16 -0.36730460
## 17 -0.50575156
## 18 -0.62074269
## -------------------------
The item ‘self_conscious’ seems to be problematic. Perhaps this item was interpreted differently by participants.
Looking at the correlation plot
corr_data = mood_efa_nodrug %>%
select(-c(participant, session, stage, task_type, Drug, Stress, task_order))
corr = round(cor(corr_data, use = "complete.obs"),2)
pmat = cor_pmat(corr_data)
ggcorrplot(corr, hc.order = T, lab = F)
I decided to use the factors yielded by oblimin, in only placebo sessions, with ‘self_conscious’ removed.
Factor based scores
positive_affect = c('happy', 'good', 'confident', 'safe', 'relaxed')
negative_affect = c('stressed', 'anxious', 'vulnerable', 'embarassed', 'distressed')
stress_activation = c('warm_face', 'heart_palpitations', 'shaking')
meh = c('bored', 'indifferent', 'dull')
# add basline correction
mood_with_bsln = mood_data_c %>%
arrange(participant, session, stage, item_name) %>%
ddply(.(participant, session, item_name), transform, bsln_corr = response - response[3])
mood_data_with_factors = mood_with_bsln %>%
mutate(factor = ifelse(item_name %in% positive_affect, 'positive_affect',
ifelse(item_name %in% negative_affect, 'negative_affect',
ifelse(item_name %in% stress_activation, 'stress_activation',
ifelse(item_name %in% meh, 'meh',
'NA'))))) %>%
filter(factor != 'NA') %>%
group_by(participant, session, stage, factor) %>%
dplyr::mutate(factor_score = mean(response), factor_score_bsln_corr = mean(bsln_corr)) %>%
ungroup()
# positive_affect = c('happy', 'good', 'confident', 'safe', 'relaxed')
# negative_affect = c('shaking', 'stressed', 'anxious', 'heart_palpitations',
# 'vulnerable', 'embarassed', 'distressed', 'warm_face')
# meh = c('bored', 'indifferent', 'dull')
#
# mood_data_with_factors = mood_data_c %>%
# mutate(factor = ifelse(item_name %in% positive_affect, 'positive_affect',
# ifelse(item_name %in% negative_affect, 'negative_affect',
# ifelse(item_name %in% meh, 'meh',
# 'NA')))) %>%
# filter(factor != 'NA') %>%
# group_by(participant, session, stage, factor) %>%
# dplyr::mutate(factor_score = mean(response))
mood_factors = mood_data_with_factors %>%
pivot_wider(id_cols = c(participant, session, stage, task_type, Stress, Drug, task_order),
names_from = factor, values_from = factor_score,
values_fn = mean)
mood_factors_bsln_corr = mood_data_with_factors %>%
pivot_wider(id_cols = c(participant, session, stage, task_type, Stress, Drug, task_order),
names_from = factor, values_from = factor_score_bsln_corr,
values_fn = mean) %>%
filter(participant %!in% c(2,29))
# mood_factors = mood_factors %>%
# arrange(participant, session, stage) %>%
# ddply(.(participant, session), transform, neg_bsln_corr = negative_affect - negative_affect[3])%>%
# ddply(.(participant, session), transform, pos_bsln_corr = positive_affect - positive_affect[3])%>%
# ddply(.(participant, session), transform, act_bsln_corr = stress_activation - stress_activation[3])%>%
# ddply(.(participant, session), transform, meh_bsln_corr = meh - meh[3])
mood_factors_placebo = mood_factors %>%
filter(Drug == 'placebo')
Baseline corrected plots
#here is the basline corrected factor data:
mood_factors_bsln_corr
#looks good
#plot the raw values:
mood_factors_bsln_corr %>%
filter(session == 1) %>%
ggplot(aes(x = stage, y = negative_affect, color = participant, group = participant)) +
geom_line()
#okay so it looks like it's doing what i want it to do
#it has actually been baseline corrected
#so the question now is how to get the means + ci from bsln corr data
# mood_factors_bsln_corr %>%
# filter(stage == 3 & negative_affect !=0)
#
# mood_factors_bsln_corr %>%
# filter(ifelse(stage == 3 & negative_affect !=0,F,T))
mood_factors_bsln_corr %>%
filter(participant != 2 & session != 2) %>%
summarySEwithin2(measurevar = 'negative_affect',
withinvars = c('task_type', 'stage'), idvar = 'participant') %>%
ggplot(aes(x = stage, y = negative_affect, color = task_type, group = c(task_type))) +
geom_line()
mood_factors_bsln_corr %>%
filter(participant != 2 & session != 2) %>%
summarySEwithin2(measurevar = 'negative_affect',
withinvars = c('task_type', 'stage'), idvar = 'participant') %>%
ggplot(aes(x = stage, y = negative_affect, color = task_type, group = c(task_type))) +
geom_line() +
geom_ribbon(aes(ymin = negative_affect-ci, ymax = negative_affect+ci, color = task_type, fill = task_type), alpha = 0.2)
Factor based score diagnostics
mood_factors_bsln_corr %>%
ggplot(aes(x=negative_affect)) +
geom_histogram() +
labs(title = "baseline corrected negative affect distribution")
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
mood_factors_bsln_corr %>%
ggplot(aes(x=positive_affect)) +
geom_histogram() +
labs(title = "baseline corrected positive affect distribution")
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
mood_factors_bsln_corr %>%
ggplot(aes(x=stress_activation)) +
geom_histogram() +
labs(title = "baseline corrected stress activation distribution")
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
model <- lm(negative_affect ~ Stress*stage, data = mood_factors_bsln_corr)
# Create a QQ plot of residuals
ggqqplot(residuals(model)) +
labs(title = "negative affect qqplot")
model <- lm(positive_affect ~ Stress*stage, data = mood_factors_bsln_corr)
# Create a QQ plot of residuals
ggqqplot(residuals(model)) +
labs(title = "positive affect qqplot")
model <- lm(stress_activation ~ Stress*stage, data = mood_factors_bsln_corr)
# Create a QQ plot of residuals
ggqqplot(residuals(model)) +
labs(title = "stress activation qqplot")
All sessions
positive_means = mood_factors%>%
summarySEwithin2(measurevar = 'positive_affect', withinvars = c('Stress', 'stage'), idvar = 'participant') %>%
mutate(lower = positive_affect-ci, upper = positive_affect+ci)
negative_means = mood_factors%>%
summarySEwithin2(measurevar = 'negative_affect', withinvars = c('Stress', 'stage'), idvar = 'participant') %>%
mutate(lower = negative_affect-ci, upper = negative_affect+ci)
activation_means = mood_factors%>%
summarySEwithin2(measurevar = 'stress_activation', withinvars = c('Stress', 'stage'), idvar = 'participant') %>%
mutate(lower = stress_activation-ci, upper = stress_activation+ci)
meh_means = mood_factors%>%
summarySEwithin2(measurevar = 'meh', withinvars = c('Stress', 'stage'), idvar = 'participant') %>%
mutate(lower = meh-ci, upper = meh+ci)
positive_means %>%
ggplot(aes(x = stage, y = positive_affect, color = Stress, group = Stress)) +
geom_line() +
geom_ribbon(aes(ymin = lower, ymax = upper, color = Stress, fill = Stress), alpha = .2) +
labs(title = "Positive Affect - all sessions") +
scale_x_discrete(labels = c("Bsl 1", "Bsl 2", "Bsl 3", "Induction", "Drug 1",
"Reminder 1", "Self-admin", "Reminder 2", "Bandit",
"Traffic", "Drug 2")) +
theme(axis.text.x = element_text(angle = 45, hjust=1), axis.text=element_text(size=12),
axis.title=element_text(size=14,face="bold"), plot.title = element_text(size=20, face = "bold"),
legend.title=element_text(size=10), legend.text=element_text(size=9))
negative_means %>%
ggplot(aes(x = stage, y = negative_affect, color = Stress, group = Stress)) +
geom_line() +
geom_ribbon(aes(ymin = lower, ymax = upper, color = Stress, fill = Stress), alpha = .2) +
labs(title = "Negative Affect - all sessions") +
scale_x_discrete(labels = c("Bsl 1", "Bsl 2", "Bsl 3", "Induction", "Drug 1",
"Reminder 1", "Self-admin", "Reminder 2", "Bandit",
"Traffic", "Drug 2")) +
theme(axis.text.x = element_text(angle = 45, hjust=1), axis.text=element_text(size=12),
axis.title=element_text(size=14,face="bold"), plot.title = element_text(size=20, face = "bold"),
legend.title=element_text(size=10), legend.text=element_text(size=9))
activation_means %>%
ggplot(aes(x = stage, y = stress_activation, color = Stress, group = Stress)) +
geom_line() +
geom_ribbon(aes(ymin = lower, ymax = upper, color = Stress, fill = Stress), alpha = .2) +
labs(title = "Stress Activation - all sessions") +
scale_x_discrete(labels = c("Bsl 1", "Bsl 2", "Bsl 3", "Induction", "Drug 1",
"Reminder 1", "Self-admin", "Reminder 2", "Bandit",
"Traffic", "Drug 2")) +
theme(axis.text.x = element_text(angle = 45, hjust=1), axis.text=element_text(size=12),
axis.title=element_text(size=14,face="bold"), plot.title = element_text(size=20, face = "bold"),
legend.title=element_text(size=10), legend.text=element_text(size=9))
meh_means %>%
ggplot(aes(x = stage, y = meh, color = Stress, group = Stress)) +
geom_line() +
geom_ribbon(aes(ymin = lower, ymax = upper, color = Stress, fill = Stress), alpha = .2) +
labs(title = "meh - all sessions") +
scale_x_discrete(labels = c("Bsl 1", "Bsl 2", "Bsl 3", "Induction", "Drug 1",
"Reminder 1", "Self-admin", "Reminder 2", "Bandit",
"Traffic", "Drug 2")) +
theme(axis.text.x = element_text(angle = 45, hjust=1), axis.text=element_text(size=12),
axis.title=element_text(size=14,face="bold"), plot.title = element_text(size=20, face = "bold"),
legend.title=element_text(size=10), legend.text=element_text(size=9))
positive_means = mood_factors%>%
summarySEwithin2(measurevar = 'positive_affect', withinvars = c('task_type', 'stage'), idvar = 'participant') %>%
mutate(lower = positive_affect-ci, upper = positive_affect+ci)
negative_means = mood_factors%>%
summarySEwithin2(measurevar = 'negative_affect', withinvars = c('task_type', 'stage'), idvar = 'participant') %>%
mutate(lower = negative_affect-ci, upper = negative_affect+ci)
activation_means = mood_factors%>%
summarySEwithin2(measurevar = 'stress_activation', withinvars = c('task_type', 'stage'), idvar = 'participant') %>%
mutate(lower = stress_activation-ci, upper = stress_activation+ci)
meh_means = mood_factors%>%
summarySEwithin2(measurevar = 'meh', withinvars = c('task_type', 'stage'), idvar = 'participant') %>%
mutate(lower = meh-ci, upper = meh+ci)
positive_means %>%
ggplot(aes(x = stage, y = positive_affect, color = task_type, group = task_type)) +
geom_line() +
geom_ribbon(aes(ymin = lower, ymax = upper, color = task_type, fill = task_type), alpha = .2) +
labs(title = "Positive Affect - all sessions") +
scale_x_discrete(labels = c("Bsl 1", "Bsl 2", "Bsl 3", "Induction", "Drug 1",
"Reminder 1", "Self-admin", "Reminder 2", "Bandit",
"Traffic", "Drug 2")) +
theme(axis.text.x = element_text(angle = 45, hjust=1), axis.text=element_text(size=12),
axis.title=element_text(size=14,face="bold"), plot.title = element_text(size=20, face = "bold"),
legend.title=element_text(size=10), legend.text=element_text(size=9))
negative_means %>%
ggplot(aes(x = stage, y = negative_affect, color = task_type, group = task_type)) +
geom_line() +
geom_ribbon(aes(ymin = lower, ymax = upper, color = task_type, fill = task_type), alpha = .2) +
labs(title = "Negative Affect - all sessions") +
scale_x_discrete(labels = c("Bsl 1", "Bsl 2", "Bsl 3", "Induction", "Drug 1",
"Reminder 1", "Self-admin", "Reminder 2", "Bandit",
"Traffic", "Drug 2")) +
theme(axis.text.x = element_text(angle = 45, hjust=1), axis.text=element_text(size=12),
axis.title=element_text(size=14,face="bold"), plot.title = element_text(size=20, face = "bold"),
legend.title=element_text(size=10), legend.text=element_text(size=9))
activation_means %>%
ggplot(aes(x = stage, y = stress_activation, color = task_type, group = task_type)) +
geom_line() +
geom_ribbon(aes(ymin = lower, ymax = upper, color = task_type, fill = task_type), alpha = .2) +
labs(title = "Stress Activation - all sessions") +
scale_x_discrete(labels = c("Bsl 1", "Bsl 2", "Bsl 3", "Induction", "Drug 1",
"Reminder 1", "Self-admin", "Reminder 2", "Bandit",
"Traffic", "Drug 2")) +
theme(axis.text.x = element_text(angle = 45, hjust=1), axis.text=element_text(size=12),
axis.title=element_text(size=14,face="bold"), plot.title = element_text(size=20, face = "bold"),
legend.title=element_text(size=10), legend.text=element_text(size=9))
meh_means %>%
ggplot(aes(x = stage, y = meh, color = task_type, group = task_type)) +
geom_line() +
geom_ribbon(aes(ymin = lower, ymax = upper, color = task_type, fill = task_type), alpha = .2) +
labs(title = "meh - all sessions") +
scale_x_discrete(labels = c("Bsl 1", "Bsl 2", "Bsl 3", "Induction", "Drug 1",
"Reminder 1", "Self-admin", "Reminder 2", "Bandit",
"Traffic", "Drug 2")) +
theme(axis.text.x = element_text(angle = 45, hjust=1), axis.text=element_text(size=12),
axis.title=element_text(size=14,face="bold"), plot.title = element_text(size=20, face = "bold"),
legend.title=element_text(size=10), legend.text=element_text(size=9))
Placebo sessions
positive_means = mood_factors_placebo %>%
summarySEwithin2(measurevar = 'positive_affect', withinvars = c('Stress', 'stage'), idvar = 'participant') %>%
mutate(lower = positive_affect-ci, upper = positive_affect+ci)
negative_means = mood_factors_placebo %>%
summarySEwithin2(measurevar = 'negative_affect', withinvars = c('Stress', 'stage'), idvar = 'participant') %>%
mutate(lower = negative_affect-ci, upper = negative_affect+ci)
activation_means = mood_factors_placebo %>%
summarySEwithin2(measurevar = 'stress_activation', withinvars = c('Stress', 'stage'), idvar = 'participant') %>%
mutate(lower = stress_activation-ci, upper = stress_activation+ci)
meh_means = mood_factors_placebo %>%
summarySEwithin2(measurevar = 'meh', withinvars = c('Stress', 'stage'), idvar = 'participant') %>%
mutate(lower = meh-ci, upper = meh+ci)
positive_means %>%
ggplot(aes(x = stage, y = positive_affect, color = Stress, group = Stress)) +
geom_line() +
geom_ribbon(aes(ymin = lower, ymax = upper, color = Stress, fill = Stress), alpha = .2) +
labs(title = "Positive Affect - placebo sessions") +
scale_x_discrete(labels = c("Bsl 1", "Bsl 2", "Bsl 3", "Induction", "Drug 1",
"Reminder 1", "Self-admin", "Reminder 2", "Bandit",
"Traffic", "Drug 2")) +
theme(axis.text.x = element_text(angle = 45, hjust=1), axis.text=element_text(size=12),
axis.title=element_text(size=14,face="bold"), plot.title = element_text(size=20, face = "bold"),
legend.title=element_text(size=10), legend.text=element_text(size=9))
negative_means %>%
ggplot(aes(x = stage, y = negative_affect, color = Stress, group = Stress)) +
geom_line() +
geom_ribbon(aes(ymin = lower, ymax = upper, color = Stress, fill = Stress), alpha = .2) +
labs(title = "Negative Affect - placebo sessions ") +
scale_x_discrete(labels = c("Bsl 1", "Bsl 2", "Bsl 3", "Induction", "Drug 1",
"Reminder 1", "Self-admin", "Reminder 2", "Bandit",
"Traffic", "Drug 2")) +
theme(axis.text.x = element_text(angle = 45, hjust=1), axis.text=element_text(size=12),
axis.title=element_text(size=14,face="bold"), plot.title = element_text(size=20, face = "bold"),
legend.title=element_text(size=10), legend.text=element_text(size=9))
activation_means %>%
ggplot(aes(x = stage, y = stress_activation, color = Stress, group = Stress)) +
geom_line() +
geom_ribbon(aes(ymin = lower, ymax = upper, color = Stress, fill = Stress), alpha = .2) +
labs(title = "Stress Activation - placebo sessions ") +
scale_x_discrete(labels = c("Bsl 1", "Bsl 2", "Bsl 3", "Induction", "Drug 1",
"Reminder 1", "Self-admin", "Reminder 2", "Bandit",
"Traffic", "Drug 2")) +
theme(axis.text.x = element_text(angle = 45, hjust=1), axis.text=element_text(size=12),
axis.title=element_text(size=14,face="bold"), plot.title = element_text(size=20, face = "bold"),
legend.title=element_text(size=10), legend.text=element_text(size=9))
meh_means %>%
ggplot(aes(x = stage, y = meh, color = Stress, group = Stress)) +
geom_line() +
geom_ribbon(aes(ymin = lower, ymax = upper, color = Stress, fill = Stress), alpha = .2) +
labs(title = "meh - placebo sessions ") +
scale_x_discrete(labels = c("Bsl 1", "Bsl 2", "Bsl 3", "Induction", "Drug 1",
"Reminder 1", "Self-admin", "Reminder 2", "Bandit",
"Traffic", "Drug 2")) +
theme(axis.text.x = element_text(angle = 45, hjust=1), axis.text=element_text(size=12),
axis.title=element_text(size=14,face="bold"), plot.title = element_text(size=20, face = "bold"),
legend.title=element_text(size=10), legend.text=element_text(size=9))
positive_means = mood_factors_placebo %>%
summarySEwithin2(measurevar = 'positive_affect', withinvars = c('task_type', 'stage'), idvar = 'participant') %>%
mutate(lower = positive_affect-ci, upper = positive_affect+ci)
negative_means = mood_factors_placebo %>%
summarySEwithin2(measurevar = 'negative_affect', withinvars = c('task_type', 'stage'), idvar = 'participant') %>%
mutate(lower = negative_affect-ci, upper = negative_affect+ci)
activation_means = mood_factors_placebo %>%
summarySEwithin2(measurevar = 'stress_activation', withinvars = c('task_type', 'stage'), idvar = 'participant') %>%
mutate(lower = stress_activation-ci, upper = stress_activation+ci)
meh_means = mood_factors_placebo %>%
summarySEwithin2(measurevar = 'meh', withinvars = c('task_type', 'stage'), idvar = 'participant') %>%
mutate(lower = meh-ci, upper = meh+ci)
positive_means %>%
ggplot(aes(x = stage, y = positive_affect, color = task_type, group = task_type)) +
geom_line() +
geom_ribbon(aes(ymin = lower, ymax = upper, color = task_type, fill = task_type), alpha = .2) +
labs(title = "Positive Affect") +
scale_color_manual(values=c("tsst"="#F8766D", "singing" = "#C77CFF", "color"="#00BFC4", "music" = "#7CAE00"))+
scale_fill_manual(values=c("tsst"="#F8766D", "singing" = "#C77CFF", "color"="#00BFC4", "music" = "#7CAE00"))+
theme_bw()+
scale_x_discrete(labels = c("Bsl 1", "Bsl 2", "Bsl 3", "Induction", "Drug 1",
"Reminder 1", "Self-admin", "Reminder 2", "Bandit",
"Traffic", "Drug 2")) +
theme(axis.text.x = element_text(angle = 45, hjust=1), axis.text=element_text(size=12),
axis.title=element_text(size=14,face="bold"), plot.title = element_text(size=20, face = "bold"),
legend.title=element_text(size=10), legend.text=element_text(size=9))
negative_means %>%
ggplot(aes(x = stage, y = negative_affect, color = task_type, group = task_type)) +
geom_line() +
geom_ribbon(aes(ymin = lower, ymax = upper, color = task_type, fill = task_type), alpha = .2) +
labs(title = "Negative Affect - placebo sessions only") +
scale_x_discrete(labels = c("Bsl 1", "Bsl 2", "Bsl 3", "Induction", "Drug 1",
"Reminder 1", "Self-admin", "Reminder 2", "Bandit",
"Traffic", "Drug 2")) +
theme(axis.text.x = element_text(angle = 45, hjust=1), axis.text=element_text(size=12),
axis.title=element_text(size=14,face="bold"), plot.title = element_text(size=20, face = "bold"),
legend.title=element_text(size=10), legend.text=element_text(size=9))
activation_means %>%
ggplot(aes(x = stage, y = stress_activation, color = task_type, group = task_type)) +
geom_line() +
geom_ribbon(aes(ymin = lower, ymax = upper, color = task_type, fill = task_type), alpha = .2) +
labs(title = "Stress Activation - placebo sessions only") +
scale_x_discrete(labels = c("Bsl 1", "Bsl 2", "Bsl 3", "Induction", "Drug 1",
"Reminder 1", "Self-admin", "Reminder 2", "Bandit",
"Traffic", "Drug 2")) +
theme(axis.text.x = element_text(angle = 45, hjust=1), axis.text=element_text(size=12),
axis.title=element_text(size=14,face="bold"), plot.title = element_text(size=20, face = "bold"),
legend.title=element_text(size=10), legend.text=element_text(size=9))
meh_means %>%
ggplot(aes(x = stage, y = meh, color = task_type, group = task_type)) +
geom_line() +
geom_ribbon(aes(ymin = lower, ymax = upper, color = task_type, fill = task_type), alpha = .2) +
labs(title = "meh - placebo sessions only") +
scale_x_discrete(labels = c("Bsl 1", "Bsl 2", "Bsl 3", "Induction", "Drug 1",
"Reminder 1", "Self-admin", "Reminder 2", "Bandit",
"Traffic", "Drug 2")) +
theme(axis.text.x = element_text(angle = 45, hjust=1), axis.text=element_text(size=12),
axis.title=element_text(size=14,face="bold"), plot.title = element_text(size=20, face = "bold"),
legend.title=element_text(size=10), legend.text=element_text(size=9))
Baseline correction plots
positive_bsl_means = mood_factors_bsln_corr %>%
summarySEwithin2(measurevar = 'positive_affect', withinvars = c('task_type', 'stage'), idvar = 'participant') %>%
mutate(lower = positive_affect-ci, upper = positive_affect+ci)
negative_bsl_means = mood_factors_bsln_corr %>%
summarySEwithin2(measurevar = 'negative_affect', withinvars = c('task_type', 'stage'), idvar = 'participant') %>%
mutate(lower = negative_affect-ci, upper = negative_affect+ci)
positive_bsl_means %>%
ggplot(aes(x = stage, y = positive_affect, color = task_type, group = task_type)) +
geom_line() +
geom_ribbon(aes(ymin = lower, ymax = upper, color = task_type, fill = task_type), alpha = .2) +
labs(title = "Positive Affect - placebo sessions, basline corrected") +
scale_x_discrete(labels = c("Bsl 1", "Bsl 2", "Bsl 3", "Induction", "Drug 1",
"Reminder 1", "Self-admin", "Reminder 2", "Bandit",
"Traffic", "Drug 2")) +
theme(axis.text.x = element_text(angle = 45, hjust=1), axis.text=element_text(size=12),
axis.title=element_text(size=14,face="bold"), plot.title = element_text(size=20, face = "bold"),
legend.title=element_text(size=10), legend.text=element_text(size=9))
negative_bsl_means %>%
ggplot(aes(x = stage, y = negative_affect, color = task_type, group = task_type)) +
geom_line() +
geom_ribbon(aes(ymin = lower, ymax = upper, color = task_type, fill = task_type), alpha = .2) +
labs(title = "Negative Affect - placebo sessions, baseline corrected") +
scale_x_discrete(labels = c("Bsl 1", "Bsl 2", "Bsl 3", "Induction", "Drug 1",
"Reminder 1", "Self-admin", "Reminder 2", "Bandit",
"Traffic", "Drug 2")) +
theme(axis.text.x = element_text(angle = 45, hjust=1), axis.text=element_text(size=12),
axis.title=element_text(size=14,face="bold"), plot.title = element_text(size=20, face = "bold"),
legend.title=element_text(size=10), legend.text=element_text(size=9))
For the first model I will only look at the state induction, i.e. Q3 & 4, for all sessions
For the second model I will look at stages 3-8 for placebo sessions only. In this model the repetition of the ReSST is a between-subjects factor (rep: regular/musical)
Select the appropriate data:
predrug_only_factors = mood_factors %>%
filter(stage %in% c(3,4))
placebo_only_factors = mood_factors_placebo %>%
filter(stage %in% c(5,6,7,8)) %>%
mutate(rep = case_when(task_type == 'color' | task_type == 'tsst' ~ 'regular',
task_type == 'music' | task_type == 'singing' ~ 'musical'))
induction_bslcorr = mood_factors_bsln_corr %>%
filter(stage %in% c(3,4),
participant %!in% c(19,29))
Some data is missing, so I’m removing those participants for anovas
predrug_only_factors_anova = predrug_only_factors %>%
filter(participant != 19 & participant != 29)
placebo_only_factors_anova = placebo_only_factors %>%
filter(participant != 19 & participant != 29)
# predrug_only_factors_anova = predrug_only_factors #%>%
# #filter(participant != 2 & participant != 58 & participant != 99)
#
# placebo_only_factors_anova = placebo_only_factors #%>%
# #filter(participant != 2 & participant != 58 & participant != 99)
Means + SE with summarySEwithin2
predrug_only_factors_anova %>%
summarySEwithin2(measurevar = 'negative_affect', withinvars = c('task_type','stage'),
idvar = 'participant')
predrug_only_factors_anova %>%
summarySEwithin2(measurevar = 'positive_affect', withinvars = c('task_type','stage'),
idvar = 'participant')
predrug_only_factors_anova %>%
summarySEwithin2(measurevar = 'stress_activation', withinvars = c('task_type','stage'),
idvar = 'participant')
predrug_only_factors_anova %>%
summarySEwithin2(measurevar = 'negative_affect', withinvars = c('Stress','stage'),
idvar = 'participant')
predrug_only_factors_anova %>%
summarySEwithin2(measurevar = 'positive_affect', withinvars = c('Stress','stage'),
idvar = 'participant')
predrug_only_factors_anova %>%
summarySEwithin2(measurevar = 'stress_activation', withinvars = c('Stress','stage'),
idvar = 'participant')
placebo_only_factors_anova %>%
summarySEwithin2(measurevar = 'negative_affect', withinvars = c('Stress','stage'),
idvar = 'participant')
placebo_only_factors_anova %>%
summarySEwithin2(measurevar = 'positive_affect', withinvars = c('Stress','stage'),
idvar = 'participant')
placebo_only_factors_anova %>%
summarySEwithin2(measurevar = 'stress_activation', withinvars = c('Stress','stage'),
idvar = 'participant')
Model 1: Induction, all sessions
Model 1a negative affect, pre-drug, all sessions
#significant interaction of task_type and stage
predrug_only_factors_anova %>%
anova_test(dv = negative_affect, wid = participant, within = c(task_type, stage), between = task_order) %>%
get_anova_table()
# outliers = predrug_only_factors_anova %>%
# identify_outliers(negative_affect) %>%
# filter(is.extreme == T)
# wo_outliers = anti_join(predrug_only_factors_anova, outliers, by = c('participant', 'session', 'stage'))
#
# wo_outliers %>%
# anova_test(dv = negative_affect, wid = participant, within = c(task_type, stage), between = task_order) %>%
# get_anova_table()
#
# wo_outliers %>%
# group_by(stage) %>%
# anova_test(dv = negative_affect, wid = participant, within = task_type) %>%
# get_anova_table() %>%
# adjust_pvalue(method = "hochberg")
#
# wo_outliers %>%
# #group_by(task_order) %>%
# anova_test(dv = negative_affect, wid = participant, within = c(task_type)) %>%
# get_anova_table() %>%
# adjust_pvalue(method = "hochberg")
#significant effect of task_type at both stages
one.way <- predrug_only_factors_anova %>%
group_by(stage) %>%
anova_test(dv = negative_affect, wid = participant, within = task_type) %>%
get_anova_table() %>%
adjust_pvalue(method = "hochberg")
one.way
#sig diff between both control and both stress tasks at stage 4, all else ns
pwc = predrug_only_factors_anova %>%
group_by(stage) %>%
pairwise_t_test(negative_affect ~ task_type, paired = TRUE,
p.adjust.method = "hochberg")
pwc
Model 1b positive affect, pre-drug, all sessions
#significant interaction of task_type and stage
predrug_only_factors_anova %>%
anova_test(dv = positive_affect, wid = participant, within = c(task_type, stage), between = task_order) %>%
get_anova_table()
# outliers = predrug_only_factors_anova %>%
# identify_outliers(positive_affect) %>%
# filter(is.extreme == T)
#
# wo_outliers = anti_join(predrug_only_factors_anova, outliers, by = c('participant', 'session', 'stage'))
#
# wo_outliers %>%
# anova_test(dv = positive_affect, wid = participant, within = c(task_type, stage), between = task_order) %>%
# get_anova_table()
#
# one.way <- wo_outliers %>%
# group_by(stage) %>%
# anova_test(dv = positive_affect, wid = participant, within = task_type) %>%
# get_anova_table() %>%
# adjust_pvalue(method = "hochberg")
# one.way
#
# pwc <- wo_outliers %>%
# group_by(stage) %>%
# pairwise_t_test(
# positive_affect ~ task_type, paired = TRUE,
# p.adjust.method = "hochberg"
# )
#
# pwc %>%
# filter(stage == 4)
#significant effect of task_type at stage 4
one.way <- predrug_only_factors_anova %>%
group_by(stage) %>%
anova_test(dv = positive_affect, wid = participant, within = task_type) %>%
get_anova_table() %>%
adjust_pvalue(method = "hochberg")
one.way
#sig diff between both control and both stress tasks at stage 4, all else ns
pwc <- predrug_only_factors_anova %>%
group_by(stage) %>%
pairwise_t_test(
positive_affect ~ task_type, paired = TRUE,
p.adjust.method = "hochberg"
)
pwc %>%
filter(stage == 4)
Model 1c stress activation, pre-drug, all sessions
#significant interaction of task_type and stage
predrug_only_factors_anova %>%
anova_test(dv = stress_activation, wid = participant, within = c(task_type, stage), between = task_order) %>%
get_anova_table()
# outliers = predrug_only_factors_anova %>%
# identify_outliers(stress_activation) %>%
# filter(is.extreme == T)
#
# wo_outliers = anti_join(predrug_only_factors_anova, outliers, by = c('participant', 'session', 'stage'))
#
# wo_outliers %>%
# anova_test(dv = stress_activation, wid = participant, within = c(task_type, stage), between = task_order) %>%
# get_anova_table()
#significant effect of task_type at stage 4
# one.way <- predrug_only_factors_anova %>%
# group_by(stage) %>%
# anova_test(dv = stress_activation, wid = participant, within = task_type) %>%
# get_anova_table() %>%
# adjust_pvalue(method = "hochberg")
# one.way
#sig diff between both control and both stress tasks at stage 4, all else ns
pwc <- predrug_only_factors_anova %>%
group_by(stage) %>%
pairwise_t_test(
stress_activation ~ task_type, paired = TRUE,
p.adjust.method = "hochberg"
)
pwc
Model 1d meh, pre-drug, all sessions
# predrug_only_factors_anova %>%
# anova_test(dv = meh, wid = participant, within = c(task_type, stage)) %>%
# get_anova_table()
#
# one.way <- predrug_only_factors_anova %>%
# group_by(task_type) %>%
# anova_test(dv = meh, wid = participant, within = stage) %>%
# get_anova_table() %>%
# adjust_pvalue(method = "bonferroni")
# one.way
#
# pwc <- predrug_only_factors_anova %>%
# group_by(task_type) %>%
# pairwise_t_test(
# meh ~ stage, paired = TRUE,
# p.adjust.method = "bonferroni"
# )
# pwc
Model 2: induction - reminder 2, placebo sessions only
Model 2a negative affect, placebo only, all stages
#significant interaction between stress and stage
placebo_only_factors_anova %>%
anova_test(dv = negative_affect, wid = participant, within = c(Stress, stage), between = c(rep, task_order)) %>%
get_anova_table()
# outliers = placebo_only_factors_anova %>%
# identify_outliers(negative_affect) %>%
# filter(is.extreme == T)
#
# wo_outliers = anti_join(placebo_only_factors_anova, outliers, by = c('participant', 'session', 'stage'))
#
# wo_outliers %>%
# anova_test(dv = negative_affect, wid = participant, within = c(Stress, stage), between = c(rep, task_order)) %>%
# get_anova_table()
# #significant effect of stress at stages 4-8
# one.way <- placebo_only_factors_anova %>%
# group_by(stage) %>%
# anova_test(dv = negative_affect, wid = participant, within = Stress) %>%
# get_anova_table() %>%
# adjust_pvalue(method = "hochberg")
# one.way
#sig diff between pleasant and stress at stages 4-8
pwc <- placebo_only_factors_anova %>%
group_by(stage) %>%
pairwise_t_test(
negative_affect ~ Stress, paired = TRUE,
p.adjust.method = "hochberg"
)
pwc %>%
filter(stage %in% c(5,6,7,8))
Model 2b positive affect, placebo only, all stages
#significant interaction between stress and stage
placebo_only_factors_anova %>%
anova_test(dv = positive_affect, wid = participant, within = c(Stress, stage), between = c(rep, task_order)) %>%
get_anova_table()
#
# outliers = placebo_only_factors_anova %>%
# identify_outliers(positive_affect) %>%
# filter(is.extreme == T)
#
# wo_outliers = anti_join(placebo_only_factors_anova, outliers, by = c('participant', 'session', 'stage'))
#
# wo_outliers %>%
# anova_test(dv = positive_affect, wid = participant, within = c(Stress, stage), between = c(rep, task_order)) %>%
# get_anova_table()
# #significant effect of Stress at stages 4-8
# one.way <- placebo_only_factors_anova %>%
# group_by(stage) %>%
# anova_test(dv = positive_affect, wid = participant, within = Stress) %>%
# get_anova_table() %>%
# adjust_pvalue(method = "hochberg")
# one.way
pwc <- placebo_only_factors_anova %>%
group_by(stage) %>%
pairwise_t_test(
positive_affect ~ Stress, paired = TRUE,
p.adjust.method = "hochberg"
)
pwc %>%
filter(stage %in% c(5,6,7,8))
#sig diff between pleasant and stress at stages 4,6,8
Model 2c Stress activation, placebo only, all stages
#significant interaction between stress and stage
placebo_only_factors_anova %>%
anova_test(dv = stress_activation, wid = participant, within = c(Stress, stage), between = c(rep, task_order)) %>%
get_anova_table()
# outliers = placebo_only_factors_anova %>%
# identify_outliers(stress_activation) %>%
# filter(is.extreme == T)
#
# wo_outliers = anti_join(placebo_only_factors_anova, outliers, by = c('participant', 'session', 'stage'))
#
# wo_outliers %>%
# anova_test(dv = stress_activation, wid = participant, within = c(Stress, stage), between = c(rep, task_order)) %>%
# get_anova_table()
#significant effect of Stress at stages 4-8
# one.way <- placebo_only_factors_anova %>%
# group_by(stage) %>%
# anova_test(dv = stress_activation, wid = participant, within = Stress) %>%
# get_anova_table() %>%
# adjust_pvalue(method = "hochberg")
# one.way
pwc <- placebo_only_factors_anova %>%
group_by(stage) %>%
pairwise_t_test(
stress_activation ~ Stress, paired = TRUE,
p.adjust.method = "hochberg"
)
pwc %>%
filter(stage %in% c(5,6,7,8))
Model 2d meh, placebo only, all stages
# placebo_only_factors_anova %>%
# anova_test(dv = meh, wid = participant, within = c(Stress, stage), between = rep) %>%
# get_anova_table()
#
# two.way <- placebo_only_factors_anova %>%
# group_by(Stress) %>%
# anova_test(dv = meh, wid = participant, within = stage, between = rep) %>%
# get_anova_table() %>%
# adjust_pvalue(method = "bonferroni")
# two.way
#
# one.way <- placebo_only_factors_anova %>%
# group_by(rep, Stress) %>%
# anova_test(dv = meh, wid = participant, within = stage) %>%
# get_anova_table() %>%
# adjust_pvalue(method = "bonferroni")
# one.way
#
# pwc <- placebo_only_factors_anova %>%
# group_by(rep, Stress) %>%
# pairwise_t_test(
# meh ~ stage, paired = TRUE,
# p.adjust.method = "bonferroni"
# )
# pwc
ANOVA models but with aov() so I can get residuals, also emmeans
#model 1a
model_1a = aov(negative_affect ~ task_type*stage*task_order+Error(participant/(task_type*stage)), data = predrug_only_factors_anova)
emmeans(model_1a,specs = ~ stage | task_type)
## Note: re-fitting model with sum-to-zero contrasts
## NOTE: Results may be misleading due to involvement in interactions
## task_type = color:
## stage emmean SE df lower.CL upper.CL
## 3 6.51 1.63 315 3.307 9.72
## 4 5.38 1.63 315 2.172 8.59
##
## task_type = music:
## stage emmean SE df lower.CL upper.CL
## 3 5.04 1.63 315 1.836 8.25
## 4 4.44 1.63 315 1.236 7.65
##
## task_type = singing:
## stage emmean SE df lower.CL upper.CL
## 3 3.95 1.63 315 0.739 7.15
## 4 33.63 1.63 315 30.419 36.83
##
## task_type = tsst:
## stage emmean SE df lower.CL upper.CL
## 3 6.23 1.63 315 3.025 9.44
## 4 36.29 1.63 315 33.081 39.50
##
## Results are averaged over the levels of: task_order
## Warning: EMMs are biased unless design is perfectly balanced
## Confidence level used: 0.95
model.pr = proj(model_1a)
predrug_only_factors_anova$neg_resi = model.pr[[3]][, "Residuals"]
predrug_only_factors_anova$neg_fit = fitted.aovlist(model_1a)
ggqqplot(predrug_only_factors_anova$neg_resi) +
labs(title = "Q-Q plot - Model 1 - Negative Affect")
predrug_only_factors_anova %>%
ggplot(aes(x=neg_fit, y=neg_resi)) +
geom_point() +
geom_hline(yintercept = 0)+
labs(title = "Residuals - Model 1 - Negative Affect", x = "Fitted", y = "Residuals") +
theme_bw()
#anova_summary(model)
#model 1b
model = aov(positive_affect ~ task_type*stage*task_order+Error(participant/(task_type*stage)), data = predrug_only_factors_anova)
emmeans(model, specs = ~ stage | task_type)
## Note: re-fitting model with sum-to-zero contrasts
## NOTE: Results may be misleading due to involvement in interactions
## task_type = color:
## stage emmean SE df lower.CL upper.CL
## 3 71.5 2.29 181 67.0 76.0
## 4 73.2 2.29 181 68.6 77.7
##
## task_type = music:
## stage emmean SE df lower.CL upper.CL
## 3 71.6 2.29 181 67.1 76.1
## 4 73.7 2.29 181 69.1 78.2
##
## task_type = singing:
## stage emmean SE df lower.CL upper.CL
## 3 74.3 2.29 181 69.7 78.8
## 4 46.8 2.29 181 42.3 51.4
##
## task_type = tsst:
## stage emmean SE df lower.CL upper.CL
## 3 71.2 2.29 181 66.7 75.8
## 4 44.6 2.29 181 40.1 49.2
##
## Results are averaged over the levels of: task_order
## Warning: EMMs are biased unless design is perfectly balanced
## Confidence level used: 0.95
model.pr = proj(model)
predrug_only_factors_anova$pos_resi = model.pr[[3]][, "Residuals"]
predrug_only_factors_anova$pos_fit = fitted.aovlist(model)
ggqqplot(predrug_only_factors_anova$pos_resi) +
labs(title = "Q-Q plot - Model 1 - Positive Affect")
predrug_only_factors_anova %>%
ggplot(aes(x=pos_fit, y=pos_resi)) +
geom_point() +
geom_hline(yintercept = 0)+
labs(title = "Residuals - Model 1 - Positive Affect", x = "Fitted", y = "Residuals") +
theme_bw()
#anova_summary(model)
#model 1c
model = aov(stress_activation ~ task_type*stage*task_order+Error(participant/(task_type*stage)), data = predrug_only_factors_anova)
emmeans(model, specs = ~ stage | task_type)
## Note: re-fitting model with sum-to-zero contrasts
## NOTE: Results may be misleading due to involvement in interactions
## task_type = color:
## stage emmean SE df lower.CL upper.CL
## 3 4.68 1.6 278 1.534 7.82
## 4 5.61 1.6 278 2.467 8.75
##
## task_type = music:
## stage emmean SE df lower.CL upper.CL
## 3 3.81 1.6 278 0.668 6.95
## 4 3.82 1.6 278 0.676 6.96
##
## task_type = singing:
## stage emmean SE df lower.CL upper.CL
## 3 3.32 1.6 278 0.175 6.46
## 4 25.69 1.6 278 22.547 28.83
##
## task_type = tsst:
## stage emmean SE df lower.CL upper.CL
## 3 5.52 1.6 278 2.381 8.66
## 4 25.84 1.6 278 22.703 28.99
##
## Results are averaged over the levels of: task_order
## Warning: EMMs are biased unless design is perfectly balanced
## Confidence level used: 0.95
model.pr = proj(model)
predrug_only_factors_anova$act_resi = model.pr[[3]][, "Residuals"]
predrug_only_factors_anova$act_fit = fitted.aovlist(model)
ggqqplot(predrug_only_factors_anova$act_resi) +
labs(title = "Q-Q plot - Model 1 - Stress Activation")
predrug_only_factors_anova %>%
ggplot(aes(x=act_fit, y=act_resi)) +
geom_point() +
geom_hline(yintercept = 0)+
labs(title = "Residuals - Model 1 - Stress Activation", x = "Fitted", y = "Residuals") +
theme_bw()
#anova_summary(model)
#model 2a
model = aov(negative_affect ~ Stress*stage*task_order*rep+Error(participant/(Stress*stage)), data = placebo_only_factors_anova)
emmeans(model, specs = ~ stage | Stress)
## Note: re-fitting model with sum-to-zero contrasts
## NOTE: Results may be misleading due to involvement in interactions
## Stress = pleasant:
## stage emmean SE df lower.CL upper.CL
## 5 5.70 1.61 192 2.53 8.87
## 6 4.21 1.61 192 1.04 7.38
## 7 4.90 1.61 192 1.73 8.07
## 8 4.27 1.61 192 1.10 7.44
##
## Stress = stress:
## stage emmean SE df lower.CL upper.CL
## 5 11.03 1.61 192 7.86 14.20
## 6 22.28 1.61 192 19.11 25.45
## 7 11.57 1.61 192 8.40 14.74
## 8 16.71 1.61 192 13.54 19.88
##
## Results are averaged over the levels of: task_order, rep
## Warning: EMMs are biased unless design is perfectly balanced
## Confidence level used: 0.95
model.pr = proj(model)
placebo_only_factors_anova$neg_resi = model.pr[[3]][, "Residuals"]
placebo_only_factors_anova$neg_fit = fitted.aovlist(model)
ggqqplot(placebo_only_factors_anova$neg_resi) +
labs(title = "Q-Q plot - Model 2 - Negative Affect")
placebo_only_factors_anova %>%
ggplot(aes(x=neg_fit, y=neg_resi)) +
geom_point() +
geom_hline(yintercept = 0)+
labs(title = "Residuals - Model 2 - Negative Affect", x = "Fitted", y = "Residuals") +
theme_bw()
#anova_summary(model)
#model 2b
model = aov(positive_affect ~ Stress*stage*task_order*rep+Error(participant/(Stress*stage)), data = placebo_only_factors_anova)
emmeans(model, specs = ~ stage | Stress)
## Note: re-fitting model with sum-to-zero contrasts
## NOTE: Results may be misleading due to involvement in interactions
## Stress = pleasant:
## stage emmean SE df lower.CL upper.CL
## 5 69.1 2.29 100 64.6 73.7
## 6 70.6 2.29 100 66.0 75.1
## 7 70.1 2.29 100 65.5 74.6
## 8 71.3 2.29 100 66.8 75.9
##
## Stress = stress:
## stage emmean SE df lower.CL upper.CL
## 5 64.5 2.29 100 60.0 69.1
## 6 57.9 2.29 100 53.3 62.4
## 7 64.7 2.29 100 60.2 69.3
## 8 61.5 2.29 100 57.0 66.1
##
## Results are averaged over the levels of: task_order, rep
## Warning: EMMs are biased unless design is perfectly balanced
## Confidence level used: 0.95
model.pr = proj(model)
placebo_only_factors_anova$pos_resi = model.pr[[3]][, "Residuals"]
placebo_only_factors_anova$pos_fit = fitted.aovlist(model)
ggqqplot(placebo_only_factors_anova$pos_resi) +
labs(title = "Q-Q plot - Model 2 - Positive Affect")
placebo_only_factors_anova %>%
ggplot(aes(x=pos_fit, y=pos_resi)) +
geom_point() +
geom_hline(yintercept = 0)+
labs(title = "Residuals - Model 2 - Positive Affect", x = "Fitted", y = "Residuals") +
theme_bw()
#anova_summary(model)
#model 2c
model = aov(stress_activation ~ Stress*stage*task_order*rep+Error(participant/(Stress*stage)), data = placebo_only_factors_anova)
emmeans(model, specs = ~ stage | Stress)
## Note: re-fitting model with sum-to-zero contrasts
## NOTE: Results may be misleading due to involvement in interactions
## Stress = pleasant:
## stage emmean SE df lower.CL upper.CL
## 5 4.68 1.28 158 2.15 7.22
## 6 3.82 1.28 158 1.29 6.36
## 7 3.82 1.28 158 1.29 6.36
## 8 3.54 1.28 158 1.00 6.08
##
## Stress = stress:
## stage emmean SE df lower.CL upper.CL
## 5 7.93 1.28 158 5.40 10.47
## 6 12.26 1.28 158 9.73 14.80
## 7 9.69 1.28 158 7.15 12.22
## 8 10.45 1.28 158 7.91 12.99
##
## Results are averaged over the levels of: task_order, rep
## Warning: EMMs are biased unless design is perfectly balanced
## Confidence level used: 0.95
model.pr = proj(model)
placebo_only_factors_anova$act_resi = model.pr[[3]][, "Residuals"]
placebo_only_factors_anova$act_fit = fitted.aovlist(model)
ggqqplot(placebo_only_factors_anova$act_resi) +
labs(title = "Q-Q plot - Model 2 - Stress Activation")
placebo_only_factors_anova %>%
ggplot(aes(x=act_fit, y=act_resi)) +
geom_point() +
geom_hline(yintercept = 0)+
labs(title = "Residuals - Model 2 - Stress Activation", x = "Fitted", y = "Residuals") +
theme_bw()
#anova_summary(model)
equi_speech = mood_factors_bsln_corr %>%
filter(stage == 4 & task_type == 'tsst')
equi_sing = mood_factors_bsln_corr %>%
filter(stage == 4 & task_type == 'singing')
setdiff(equi_speech$participant, equi_sing$participant)
## character(0)
setdiff(equi_sing$participant, equi_speech$participant)
## character(0)
library(pastecs)
##
## Attaching package: 'pastecs'
## The following objects are masked from 'package:data.table':
##
## first, last
## The following objects are masked from 'package:dplyr':
##
## first, last
## The following object is masked from 'package:tidyr':
##
## extract
stat.desc(equi_speech$negative_affect)
## nbr.val nbr.null nbr.na min max range
## 60.0000000 2.0000000 0.0000000 -17.1710527 89.2324561 106.4035088
## sum median mean SE.mean CI.mean.0.95 var
## 1781.7653478 28.7171055 29.6960891 3.1175664 6.2382360 583.1532228
## std.dev coef.var
## 24.1485656 0.8131901
stat.desc(equi_sing$negative_affect)
## nbr.val nbr.null nbr.na min max range
## 60.0000000 1.0000000 0.0000000 -2.3903501 76.2938595 78.6842096
## sum median mean SE.mean CI.mean.0.95 var
## 1749.1885975 26.6721492 29.1531433 2.7300340 5.4627854 447.1851363
## std.dev coef.var
## 21.1467524 0.7253678
cor(equi_speech$negative_affect, equi_sing$negative_affect)
## [1] 0.7014116
TOSTpaired.raw(n = 60, m1 = 29.6960891, m2 = 29.1531433, sd1 = 24.1485656,
r12 = cor(equi_speech$negative_affect, equi_sing$negative_affect),
sd2 = 21.1467524, low_eqbound = -5, high_eqbound = 5,
plot = T, verbose = T)
## Note: this function is defunct. Please use tsum_TOST instead
## TOST results:
## t-value lower bound: 2.42 p-value lower bound: 0.009
## t-value upper bound: -1.95 p-value upper bound: 0.028
## degrees of freedom : 59
##
## Equivalence bounds (raw scores):
## low eqbound: -5
## high eqbound: 5
##
## TOST confidence interval:
## lower bound 90% CI: -3.28
## upper bound 90% CI: 4.366
##
## NHST confidence interval:
## lower bound 95% CI: -4.034
## upper bound 95% CI: 5.12
##
## Equivalence Test Result:
## The equivalence test was significant, t(59) = -1.948, p = 0.0281, given equivalence bounds of -5.000 and 5.000 (on a raw scale) and an alpha of 0.05.
##
##
## Null Hypothesis Test Result:
## The null hypothesis test was non-significant, t(59) = 0.237, p = 0.813, given an alpha of 0.05.
##
## NHST: don't reject null significance hypothesis that the effect is equal to 0
## TOST: reject null equivalence hypothesis
tsum_TOST(m1 = 29.6960891, m2 = 29.1531433, sd1 = 24.1485656, n1 = 60,
sd2 = 21.1467524, n2 = 60, r12 = cor(equi_speech$negative_affect, equi_sing$negative_affect),
hypothesis = 'EQU', paired = T, low_eqbound = -5, high_eqbound = 5, eqbound_type = 'raw')
##
## Paired t-test
## Hypothesis Tested: Equivalence
## Equivalence Bounds (raw):-5.000 & 5.000
## Alpha Level:0.05
## The equivalence test was significant, t(59) = -1.948, p = 2.81e-02
## The null hypothesis test was non-significant, t(59) = 0.237, p = 8.13e-01
## NHST: don't reject null significance hypothesis that the effect is equal to zero
## TOST: reject null equivalence hypothesis
##
## TOST Results
## t SE df p.value
## t-test 0.2373503 2.287529 59 0.81320729
## TOST Lower 2.4231146 2.287529 59 0.00923763
## TOST Upper -1.9484140 2.287529 59 0.02806331
##
## Effect Sizes
## estimate SE lower.ci upper.ci conf.level
## Raw 0.54294580 2.28752932 -3.2797285 4.3656201 0.9
## Hedges' g(z) -0.03044655 0.07777254 -0.1962937 0.1345649 0.9
##
## Note: SMD confidence intervals are an approximation. See vignette("SMD_calcs").
Model 1: Pre-drug only, all sessions
Model 1a negative affect
# na_model_1 = lmer(negative_affect ~ Stress*stage + (1|participant) + (Stress|task_type), data = predrug_only_factors)
#
# summary(na_model_1)
Model 1b positive affect
# #pa_model_1 = lmer(positive_affect ~ Stress*stage + (1|participant) + (Stress|task_type), data = predrug_only_factors)
# #singularity error
#
# pa_model_1 = lmer(positive_affect ~ Stress*stage + (1|participant), data = predrug_only_factors)
#
# #pa_model_1 = lmer(positive_affect ~ task_type*stage + (1|participant), data = predrug_only_factors)
#
# summary(pa_model_1)
Model 2: all stages, placebo sessions only
Model 2a negative affect
# #na_model_2 = lmer(negative_affect ~ Stress*stage + (1|participant) + (Stress|task_type), data = placebo_only_factors)
# #singularity error
#
# na_model_2 = lmer(negative_affect ~ Stress*stage + (1|participant), data = placebo_only_factors)
#
# summary(na_model_2)
Model 2b positive affect
# #pa_model_2 = lmer(positive_affect ~ Stress*stage + (1|participant) + (Stress|task_type), data = placebo_only_factors)
# #singularity error
#
# pa_model_2 = lmer(positive_affect ~ Stress*stage + (1|participant), data = placebo_only_factors)
#
# summary(pa_model_2)
Nice plots Placebo sessions
mood_plots = mood_factors %>%
filter(participant != 19 & participant != 29) %>%
filter(stage %in% c(2,3,4,5,6,7,8,9)) %>%
filter(ifelse(Drug == 'placebo' | (Drug == 'oxycodone' & stage %in% c(2,3,4)),T,F))
positive_means = mood_plots %>%
summarySEwithin2(measurevar = 'positive_affect', withinvars = c('Stress', 'stage'), idvar = 'participant') %>%
mutate(lower = positive_affect-ci, upper = positive_affect+ci)
#write.csv(positive_means,file = "summaries/pos_stress_summ.csv")
negative_means = mood_plots %>%
summarySEwithin2(measurevar = 'negative_affect', withinvars = c('Stress', 'stage'), idvar = 'participant') %>%
mutate(lower = negative_affect-ci, upper = negative_affect+ci)
#write.csv(negative_means,file = "summaries/neg_stress_summ.csv")
activation_means = mood_plots %>%
summarySEwithin2(measurevar = 'stress_activation', withinvars = c('Stress', 'stage'), idvar = 'participant') %>%
mutate(lower = stress_activation-ci, upper = stress_activation+ci)
#write.csv(activation_means,file = "summaries/act_stress_summ.csv")
meh_means = mood_plots %>%
summarySEwithin2(measurevar = 'meh', withinvars = c('Stress', 'stage'), idvar = 'participant') %>%
mutate(lower = meh-ci, upper = meh+ci)
positive_means %>%
ggplot(aes(x = stage, y = positive_affect, color = Stress, group = Stress)) +
geom_line() +
geom_ribbon(aes(ymin = lower, ymax = upper, color = Stress, fill = Stress), alpha = .2) +
scale_color_manual(name ="", labels=c("stress" = "Stress", "pleasant" = "Control"), values=c("stress"="#F8766D", "pleasant"="#00BFC4"))+
scale_fill_manual(name ="", labels=c("stress" = "Stress", "pleasant" = "Control"), values=c("stress"="#F8766D", "pleasant"="#00BFC4"))+
theme_bw()+
labs(title = "Positive Affect", x = "", y = "Positive Affect VAS") +
scale_x_discrete(labels = c("Baseline 1", "Baseline 2", "Post-induction",
"Pre-reminder 1", "Post-reminder 1", "Pre-reminder 2", "Post-reminder 2", "Post-unrelated task")) +
theme(axis.text.x = element_text(angle = 45, hjust=1), axis.text=element_text(size=10),
axis.title=element_text(size=12,face="bold"), plot.title = element_text(size=16, face = "bold"),
legend.title=element_text(size=10), legend.text=element_text(size=9)) +
ylim(0,100)
negative_means %>%
ggplot(aes(x = stage, y = negative_affect, color = Stress, group = Stress)) +
geom_line() +
geom_ribbon(aes(ymin = lower, ymax = upper, color = Stress, fill = Stress), alpha = .2) +
scale_color_manual(name ="", labels=c("stress" = "Stress", "pleasant" = "Control"), values=c("stress"="#F8766D", "pleasant"="#00BFC4"))+
scale_fill_manual(name ="", labels=c("stress" = "Stress", "pleasant" = "Control"), values=c("stress"="#F8766D", "pleasant"="#00BFC4"))+
theme_bw()+
labs(title = "Negative Affect", x = "", y = "Negative Affect VAS") +
scale_x_discrete(labels = c("Baseline 1", "Baseline 2", "Post-induction",
"Pre-reminder 1", "Post-reminder 1", "Pre-reminder 2", "Post-reminder 2", "Post-unrelated task"))+
theme(axis.text.x = element_text(angle = 45, hjust=1), axis.text=element_text(size=10),
axis.title=element_text(size=12,face="bold"), plot.title = element_text(size=16, face = "bold"),
legend.title=element_text(size=10), legend.text=element_text(size=9)) +
ylim(0,100)
activation_means %>%
ggplot(aes(x = stage, y = stress_activation, color = Stress, group = Stress)) +
geom_line() +
geom_ribbon(aes(ymin = lower, ymax = upper, color = Stress, fill = Stress), alpha = .2) +
scale_color_manual(name ="", labels=c("stress" = "Stress", "pleasant" = "Control"), values=c("stress"="#F8766D", "pleasant"="#00BFC4"))+
scale_fill_manual(name ="", labels=c("stress" = "Stress", "pleasant" = "Control"), values=c("stress"="#F8766D", "pleasant"="#00BFC4"))+
theme_bw()+
labs(title = "Stress Activation", x = "", y = "Stress Activation VAS") +
scale_x_discrete(labels = c("Baseline 1", "Baseline 2", "Post-induction",
"Pre-reminder 1", "Post-reminder 1", "Pre-reminder 2", "Post-reminder 2", "Post-unrelated task"))+
theme(axis.text.x = element_text(angle = 45, hjust=1), axis.text=element_text(size=10),
axis.title=element_text(size=12,face="bold"), plot.title = element_text(size=16, face = "bold"),
legend.title=element_text(size=10), legend.text=element_text(size=9)) +
ylim(0,100)
meh_means %>%
ggplot(aes(x = stage, y = meh, color = Stress, group = Stress)) +
geom_line() +
geom_ribbon(aes(ymin = lower, ymax = upper, color = Stress, fill = Stress), alpha = .2) +
scale_color_manual(name ="", labels=c("stress" = "Stress", "pleasant" = "Control"), values=c("stress"="#F8766D", "pleasant"="#00BFC4"))+
scale_fill_manual(name ="", labels=c("stress" = "Stress", "pleasant" = "Control"), values=c("stress"="#F8766D", "pleasant"="#00BFC4"))+
theme_bw()+
labs(title = "Tiredness", x = "", y = "Tiredness VAS") +
scale_x_discrete(labels = c("Baseline 1", "Baseline 2", "Post-induction",
"Pre-reminder 1", "Post-reminder 1", "Pre-reminder 2", "Post-reminder 2", "Post-unrelated task"))+
theme(axis.text.x = element_text(angle = 45, hjust=1), axis.text=element_text(size=10),
axis.title=element_text(size=12,face="bold"), plot.title = element_text(size=16, face = "bold"),
legend.title=element_text(size=10), legend.text=element_text(size=9)) +
ylim(0,100)
positive_means = mood_plots %>%
summarySEwithin2(measurevar = 'positive_affect', withinvars = c('task_type', 'stage'), idvar = 'participant') %>%
mutate(lower = positive_affect-ci, upper = positive_affect+ci)
#write.csv(positive_means,file = "summaries/pos_task_summ.csv")
negative_means = mood_plots %>%
summarySEwithin2(measurevar = 'negative_affect', withinvars = c('task_type', 'stage'), idvar = 'participant') %>%
mutate(lower = negative_affect-ci, upper = negative_affect+ci)
#write.csv(negative_means,file = "summaries/neg_task_summ.csv")
activation_means = mood_plots %>%
summarySEwithin2(measurevar = 'stress_activation', withinvars = c('task_type', 'stage'), idvar = 'participant') %>%
mutate(lower = stress_activation-ci, upper = stress_activation+ci)
#write.csv(activation_means,file = "summaries/act_task_summ.csv")
meh_means = mood_plots %>%
summarySEwithin2(measurevar = 'meh', withinvars = c('task_type', 'stage'), idvar = 'participant') %>%
mutate(lower = meh-ci, upper = meh+ci)
positive_means %>%
ggplot(aes(x = stage, y = positive_affect, color = task_type, group = task_type)) +
geom_line() +
geom_ribbon(aes(ymin = lower, ymax = upper, color = task_type, fill = task_type), alpha = .2) +
scale_color_manual(name = "", labels=c("tsst" = "Speech", "color" = "Color Control", "singing" = "Singing", "music" = "Music Control"),
values=c("tsst"="#F8766D", "singing" = "#C77CFF", "color"="#00BFC4", "music" = "#7CAE00"))+
scale_fill_manual(name = "", labels=c("tsst" = "Speech", "color" = "Color Control", "singing" = "Singing", "music" = "Music Control"),
values=c("tsst"="#F8766D", "singing" = "#C77CFF", "color"="#00BFC4", "music" = "#7CAE00"))+
theme_bw()+
labs(title = "Positive Affect", x = "", y = "Positive Affect VAS") +
scale_x_discrete(labels = c("Baseline 1", "Baseline 2", "Post-induction",
"Pre-reminder 1", "Post-reminder 1", "Pre-reminder 2", "Post-reminder 2", "Post-unrelated task"))+
theme(axis.text.x = element_text(angle = 45, hjust=1), axis.text=element_text(size=10),
axis.title=element_text(size=12,face="bold"), plot.title = element_text(size=16, face = "bold"),
legend.title=element_text(size=10), legend.text=element_text(size=9)) +
ylim(0,100)
negative_means %>%
ggplot(aes(x = stage, y = negative_affect, color = task_type, group = task_type)) +
geom_line() +
geom_ribbon(aes(ymin = lower, ymax = upper, color = task_type, fill = task_type), alpha = .2) +
scale_color_manual(name = "", labels=c("tsst" = "Speech", "color" = "Color Control", "singing" = "Singing", "music" = "Music Control"),
values=c("tsst"="#F8766D", "singing" = "#C77CFF", "color"="#00BFC4", "music" = "#7CAE00"))+
scale_fill_manual(name = "", labels=c("tsst" = "Speech", "color" = "Color Control", "singing" = "Singing", "music" = "Music Control"),
values=c("tsst"="#F8766D", "singing" = "#C77CFF", "color"="#00BFC4", "music" = "#7CAE00"))+
theme_bw()+
labs(title = "Negative Affect", x = "", y = "Negative Affect VAS") +
scale_x_discrete(labels = c("Baseline 1", "Baseline 2", "Post-induction",
"Pre-reminder 1", "Post-reminder 1", "Pre-reminder 2", "Post-reminder 2", "Post-unrelated task"))+
theme(axis.text.x = element_text(angle = 45, hjust=1), axis.text=element_text(size=10),
axis.title=element_text(size=12,face="bold"), plot.title = element_text(size=16, face = "bold"),
legend.title=element_text(size=10), legend.text=element_text(size=9)) +
ylim(0,100)
activation_means %>%
ggplot(aes(x = stage, y = stress_activation, color = task_type, group = task_type)) +
geom_line() +
geom_ribbon(aes(ymin = lower, ymax = upper, color = task_type, fill = task_type), alpha = .2) +
scale_color_manual(name = "", labels=c("tsst" = "Speech", "color" = "Color Control", "singing" = "Singing", "music" = "Music Control"),
values=c("tsst"="#F8766D", "singing" = "#C77CFF", "color"="#00BFC4", "music" = "#7CAE00"))+
scale_fill_manual(name = "", labels=c("tsst" = "Speech", "color" = "Color Control", "singing" = "Singing", "music" = "Music Control"),
values=c("tsst"="#F8766D", "singing" = "#C77CFF", "color"="#00BFC4", "music" = "#7CAE00"))+
theme_bw()+
labs(title = "Stress Activation", x = "", y = "Stress Activation VAS") +
scale_x_discrete(labels = c("Baseline 1", "Baseline 2", "Post-induction",
"Pre-reminder 1", "Post-reminder 1", "Pre-reminder 2", "Post-reminder 2", "Post-unrelated task"))+
theme(axis.text.x = element_text(angle = 45, hjust=1), axis.text=element_text(size=10),
axis.title=element_text(size=12,face="bold"), plot.title = element_text(size=16, face = "bold"),
legend.title=element_text(size=10), legend.text=element_text(size=9)) +
ylim(0,100)
meh_means %>%
ggplot(aes(x = stage, y = meh, color = task_type, group = task_type)) +
geom_line() +
geom_ribbon(aes(ymin = lower, ymax = upper, color = task_type, fill = task_type), alpha = .2) +
scale_color_manual(name = "", labels=c("tsst" = "Speech", "color" = "Color Control", "singing" = "Singing", "music" = "Music Control"),
values=c("tsst"="#F8766D", "singing" = "#C77CFF", "color"="#00BFC4", "music" = "#7CAE00"))+
scale_fill_manual(name = "", labels=c("tsst" = "Speech", "color" = "Color Control", "singing" = "Singing", "music" = "Music Control"),
values=c("tsst"="#F8766D", "singing" = "#C77CFF", "color"="#00BFC4", "music" = "#7CAE00"))+
theme_bw()+
labs(title = "Tiredness", x = "", y = "Tiredness VAS") +
scale_x_discrete(labels = c("Baseline 1", "Baseline 2", "Post-induction",
"Pre-reminder 1", "Post-reminder 1", "Pre-reminder 2", "Post-reminder 2", "Post-unrelated task"))+
theme(axis.text.x = element_text(angle = 45, hjust=1), axis.text=element_text(size=10),
axis.title=element_text(size=12,face="bold"), plot.title = element_text(size=16, face = "bold"),
legend.title=element_text(size=10), legend.text=element_text(size=9)) +
ylim(0,100)
okay, these are the baseline corrected plots that work:
blsn_corr_plots = mood_factors_bsln_corr %>%
filter(participant != 19 & participant != 29) %>%
filter(stage %in% c(2,3,4,5,6,7,8,9)) %>%
filter(ifelse(Drug == 'placebo' | (Drug == 'oxycodone' & stage %in% c(2,3,4)),T,F))
#Neg affect
blsn_corr_plots %>%
summarySEwithin2(measurevar = 'negative_affect', withinvars = c('task_type', 'stage'), idvar = 'participant') %>%
ggplot(aes(x = stage, y = negative_affect, color = task_type, group = task_type)) +
geom_line() +
geom_ribbon(aes(ymin = negative_affect-ci, ymax = negative_affect+ci, color = task_type, fill = task_type), alpha = .2) +
scale_color_manual(name = "", labels=c("tsst" = "Speech", "color" = "Color Control", "singing" = "Singing", "music" = "Music Control"),
values=c("tsst"="#F8766D", "singing" = "#C77CFF", "color"="#00BFC4", "music" = "#7CAE00"))+
scale_fill_manual(name = "", labels=c("tsst" = "Speech", "color" = "Color Control", "singing" = "Singing", "music" = "Music Control"),
values=c("tsst"="#F8766D", "singing" = "#C77CFF", "color"="#00BFC4", "music" = "#7CAE00")) +
theme_bw() +
labs(title = "Baseline-corrected Negative Affect", x = "", y = "Negative Affect VAS\n(baseline-corrected)")+
scale_x_discrete(labels = c("Baseline 1", "Baseline 2", "Post-induction",
"Pre-reminder 1", "Post-reminder 1", "Pre-reminder 2", "Post-reminder 2", "Post-unrelated task"))+
theme(axis.text.x = element_text(angle = 45, hjust=1), axis.text=element_text(size=10),
axis.title=element_text(size=12,face="bold"), plot.title = element_text(size=16, face = "bold"),
legend.title=element_text(size=10), legend.text=element_text(size=9))+
coord_cartesian(ylim = c(-50, 50))
#Pos affect
blsn_corr_plots %>%
summarySEwithin2(measurevar = 'positive_affect', withinvars = c('task_type', 'stage'), idvar = 'participant') %>%
ggplot(aes(x = stage, y = positive_affect, color = task_type, group = task_type)) +
geom_line() +
geom_ribbon(aes(ymin = positive_affect-ci, ymax = positive_affect+ci, color = task_type, fill = task_type), alpha = .2) +
scale_color_manual(name = "", labels=c("tsst" = "Speech", "color" = "Color Control", "singing" = "Singing", "music" = "Music Control"),
values=c("tsst"="#F8766D", "singing" = "#C77CFF", "color"="#00BFC4", "music" = "#7CAE00"))+
scale_fill_manual(name = "", labels=c("tsst" = "Speech", "color" = "Color Control", "singing" = "Singing", "music" = "Music Control"),
values=c("tsst"="#F8766D", "singing" = "#C77CFF", "color"="#00BFC4", "music" = "#7CAE00")) +
theme_bw() +
labs(title = "Baseline-corrected Positive Affect", x = "", y = "Positive Affect VAS\n(baseline-corrected)")+
scale_x_discrete(labels = c("Baseline 1", "Baseline 2", "Post-induction",
"Pre-reminder 1", "Post-reminder 1", "Pre-reminder 2", "Post-reminder 2", "Post-unrelated task"))+
theme(axis.text.x = element_text(angle = 45, hjust=1), axis.text=element_text(size=10),
axis.title=element_text(size=12,face="bold"), plot.title = element_text(size=16, face = "bold"),
legend.title=element_text(size=10), legend.text=element_text(size=9))+
coord_cartesian(ylim = c(-50, 50))
#Stress activation
blsn_corr_plots %>%
summarySEwithin2(measurevar = 'stress_activation', withinvars = c('task_type', 'stage'), idvar = 'participant') %>%
ggplot(aes(x = stage, y = stress_activation, color = task_type, group = task_type)) +
geom_line() +
geom_ribbon(aes(ymin = stress_activation-ci, ymax = stress_activation+ci, color = task_type, fill = task_type), alpha = .2) +
scale_color_manual(name = "", labels=c("tsst" = "Speech", "color" = "Color Control", "singing" = "Singing", "music" = "Music Control"),
values=c("tsst"="#F8766D", "singing" = "#C77CFF", "color"="#00BFC4", "music" = "#7CAE00"))+
scale_fill_manual(name = "", labels=c("tsst" = "Speech", "color" = "Color Control", "singing" = "Singing", "music" = "Music Control"),
values=c("tsst"="#F8766D", "singing" = "#C77CFF", "color"="#00BFC4", "music" = "#7CAE00")) +
theme_bw() +
labs(title = "Baseline-corrected Stress Activation", x = "", y = "Stress Activation VAS\n(baseline-corrected)")+
scale_x_discrete(labels = c("Baseline 1", "Baseline 2", "Post-induction",
"Pre-reminder 1", "Post-reminder 1", "Pre-reminder 2", "Post-reminder 2", "Post-unrelated task"))+
theme(axis.text.x = element_text(angle = 45, hjust=1), axis.text=element_text(size=10),
axis.title=element_text(size=12,face="bold"), plot.title = element_text(size=16, face = "bold"),
legend.title=element_text(size=10), legend.text=element_text(size=9))+
coord_cartesian(ylim = c(-50, 50))
#Tiredeness
blsn_corr_plots %>%
summarySEwithin2(measurevar = 'meh', withinvars = c('task_type', 'stage'), idvar = 'participant') %>%
ggplot(aes(x = stage, y = meh, color = task_type, group = task_type)) +
geom_line() +
geom_ribbon(aes(ymin = meh-ci, ymax = meh+ci, color = task_type, fill = task_type), alpha = .2) +
scale_color_manual(name = "", labels=c("tsst" = "Speech", "color" = "Color Control", "singing" = "Singing", "music" = "Music Control"),
values=c("tsst"="#F8766D", "singing" = "#C77CFF", "color"="#00BFC4", "music" = "#7CAE00"))+
scale_fill_manual(name = "", labels=c("tsst" = "Speech", "color" = "Color Control", "singing" = "Singing", "music" = "Music Control"),
values=c("tsst"="#F8766D", "singing" = "#C77CFF", "color"="#00BFC4", "music" = "#7CAE00")) +
theme_bw() +
labs(title = "Baseline-corrected Tiredness", x = "", y = "Tiredness VAS\n(baseline-corrected)")+
scale_x_discrete(labels = c("Baseline 1", "Baseline 2", "Post-induction",
"Pre-reminder 1", "Post-reminder 1", "Pre-reminder 2", "Post-reminder 2", "Post-unrelated task"))+
theme(axis.text.x = element_text(angle = 45, hjust=1), axis.text=element_text(size=10),
axis.title=element_text(size=12,face="bold"), plot.title = element_text(size=16, face = "bold"),
legend.title=element_text(size=10), legend.text=element_text(size=9))+
coord_cartesian(ylim = c(-50, 50))